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Abstract 

We show the impact of the electroweak, and in one instance the QCD, one-loop corrections 
on the rehc density of dark matter in the MSSM which is provided by the hghtest neutrahno. 
We cover here some of the most important scenarii: annihilation into fermions for a bino-like 
neutralino, annihilation involving gauge bosons in the case of a mixed neutralino, the neutralino- 
stau co-annihilation region and annihilation into a bottom quark pair. The corrections can 
be large and should be taken into account in view of the present and forthcoming increasing 
precision on the relic density measurements. Our calculations are made possible thanks to a 
newly developed automatic tool for the calculation at one-loop of any process in the MSSM. 
We have implemented a complete on-shell gauge invariant renormalisation scheme, with the 
possibility of switching to other schemes. In particular we will report on the impact of different 
renormalisation schemes for tan /3. 
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1 Introduction 



The last few years have witnessed spectacular advances in cosmology and astrophysics 
confirming with an unprecedented level of accuracy that ordinary matter is a minute part 
of what constitutes the Universe at large. At the same time as the LHC will be gathering 
data, a host of non coUider experiments will be carried out in search of Dark Matter, 
DM, (either direct or indirect) as well as through ever more precise determination of the 
cosmological parameters. In this new paradigm, the search for DM at the LHC is high on 
the agenda, as is of course the search for the Higgs. In fact these may be two facets of the 
New Physics that provides a resolution to the hierarchy problem posed by the Higgs in the 
Standard Model, SM.. The epitome of this New Physics is supersymmetry which among 
many advantages furnishes a good DM candidate through the lightest neutrahno, Xi- K 
future colliders discover supcrsymmetric particles and probe their properties, one could 
predict the dark matter density of the Universe and would constrain cosmology with the 
help of precision data [1,2] provided by WMAP [3] and PLANCK [4]. It would be highly 
exciting if the precision reconstruction of the relic density from observables at the colliders 
does not match PLANCK's determination, this would mean that the post-inflation era 
is most probably not entirely radiation dominated [5]. Already now the accuracy on the 
relic dark matter density is about 10% from WMAP and will soon be improved to about 
2% from PLANCK. Such level of accuracy must be matched by precise theoretical calcu- 
lations. Prom the particle physics point of view this means precision calculations of the 
annihilation and co- annihilation cross sections at least at one-loop. Quite sophisticated 
codes now exist [6,7] for the calculation of the relic density, however they arc essentially 
based on tree-level cross sections with the inclusion of some higher order effects essen- 
tially through some running couplings, masses or some effective couplings. Some of these 
corrections [6] have already been shown to be essential like the corrections to the Higgs 
couplings that can completely change the picture in the so-called Higgs funnel (annihila- 
tion mainly through the pseudo-scalar Higgs, A). The use of other approximations needs 
to be justified by complete higher order calculations which contain more than just the 
efi'ect of efi'ective couplings. In a word, the level of accuracy that will soon be reached 
requires that one is ready to tackle in a general way a full calculation at one-loop for any 
annihilation (or co-annihilation) of the neutralinos in supersymmetry, just as what one 
has been doing for the cross sections at the colliders. 

The aim of this letter is to report on the progress towards automatisation of these calcu- 
lations and to show and discuss some results on the one-loop corrected annihilation and 
co-annihilation cross sections of the LSP neutralino in the MSSM. In particular, we study 
here some of the most important scenarii: i) annihilation in the so-called bulk region 
into fermions for a bino-like neutralino, ii) co-annihilation involving the neutralino and 
the lightest stau fi, iii) annihilation into a pair of massive gauge bosons in the case of a 
mixed neutralino and iv) annihilation into bb where the pseudo-scalar Higgs pole can play 
a role. We concentrate on the electroweak corrections, although iv) is an excuse to show 
how we handle some classes of QCD corrections. 

The couple of very recent calculations of loop corrections to the relic density tackled ei- 
ther QCD corrections in extreme, though highly interesting, scenarii such as annihilation 
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into top quarks at threshold [8] and the nice study of stop-neutrahno co-annihilation [9]j. 
Some important non-perturbative electroweak effects of the Coulomb-Sommerfeld type 
that occur for TeV winos or higgsinos with a relative mass splitting between the lightest 
super symmetric particle (LSP) and the next-to-lightest supersymmetric particle (NLSP) 
of O(10~^) have been reported in [15,16]. Let us also note that, though not to be seen 
as radiative corrections to the annihilation cross sections, the temperature corrections to 
the relic density have been considered and found to be totally negligible at the level of 
10~^ [17]. A better simulation of the cosmological equation of state to derive the effec- 
tive number of relativistic degrees of freedom has been done giving corrections ranging 
from 1.5% to 3.5% [18] compared to the usual treatment as done in DarkSUSY [7] or 
micrOMEGAs [6]. 



2 General set-up and details of the calculation 

2.1 Set-up of the automatic calculation of the cross sections 

Even in the SAi, one-loop calculations of 2 ^ 2 processes involve hundreds of diagrams 
and a hand calculation is practically impracticable. Efficient automatic codes for any 
generic 2^2 process, that have now been exploited for many 2 — 3 [19,20] and even 
some 2—^4 [21,22] processes, are almost unavoidable for such calculations. For the 
electroweak theory these are the GRACE-loop [23] code and the bundle of packages based 
on FeynArts [24], FormCalc [25] and LoopTools [26], that we will refer to as FFL for 
short. 

With its much larger particle content, far greater number of parameters and more complex 
structure, the need for an automatic code at one-loop for the minimal supersymmetric 
standard model is even more of a must. A few parts that are needed for such a code have 
been developed based on an extension of [27] but, as far as we know, no complete code 
exists or is, at least publicly, available. Grace-susy [28] is now also being developed at 
one-loop and many results exist [29]. One of the main difficulties that has to be tackled 
is the implementation of the model file, since this requires that one enters the thousands 
of vertices that define the Feynman rules. On the theory side a proper renormalisation 
scheme needs to be set up, which then means extending many of these rules to include 
counterterms. When this is done one can just use, or hope to use, the machinery developed 
for the SAi, in particular the symbolic manipulation part and most importantly the loop 
integral routines including tensor reduction algorithms or any other efficient set of basis 
integrals. 

The results we will report are based on the development of a new automatic tool that 
uses and adapts modules, many of which, but not all, are part of other codes like FFL. 
This is the package SloopS whose main components and architecture we briefly sketch. 

In this application we combine LANHEP [30] (originally part of the package COMPHEP [31]) 
with the FFL bundle but with an extended and adapted LoopTools [11]. LANHEP is a 
very powerful routine that automatically generates all the sets of Feynman rules of a 

^We do not list here loop induced annihilation processes such XiXi ~^ 99 [10-13]. A very recent 
paper discusses the QCD correction to annihilation into bb in the funnel [14], however the bulk of all 
contributions has been known for sometime and implemented in micrOMEGAs already. 
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given model, the latter being defined in a simple and compact format very similar to the 
canonical coordinate representation. Use of multiplets and the super potential is built- 
in to minimize human error. The ghost Lagrangian is derived directly from the BRST 
transformations. The LANHEP module also allows to shift fields and parameters and thus 
generates counterterms most efficiently. Understandably the LANHEP output file must be 
in the format of the model file of the code it is interfaced with. In the case of FeynArts 
both the generic (Lorentz structure) and classes (particle content) files had to be given. 
Moreover because we use a non-linear gauge fixing condition [23], see below, the FeynArts 
default generic file had to be extended. 

2.2 Renormalisation and renormalisation schemes 

In the last half decade there has been an upsurge and fiurry of activity constraining mod- 
els of supersymmetry and other New Physics with the limit on the relic density delimiting 
most of the parameter space of these models. All these investigations are based on tree- 
level, sometimes with improved effective couplings, estimates of the relic density. Only in 
the last few months have some investigations [32], within mSUGRA, added a theoretical 
error estimate of C(10%), i.e, of the same order as the current experimental error. In 
these analyses based on renormalisation group running, a substantial uncertainty is due 
to the impact of the running itself [1,33]. Even if the weak scale spectrum is known, loop 
corrections to the cross sections are needed. In fact the precision one-loop calculations we 
are carrying will be most useful when confronting a measurement of the relic density once 
the microscopic properties of dark matter would have been pinned down at the collider 
and in direct/indirect detection. Henceforth we rely on the physical masses of the SUSY 
particles with the addition of some physical observables to fully reconstruct the model. 
We therefore work, as far as possible, within an on-shell scheme generalising what is done 
for the electroweak standard model. 

i) The Standard Model parameters: the fermion masses as well as the mass of the W and 
the Z are taken as input physical masses. The electric charge is defined in the Thom- 
son limit, see for example [23]. Because we are calculating corrections to processes at a 
scale 2m^o ~ 2Mz, the effect from the running electromagnetic coupling due to the light 
fermion masses will, alone, rescale the tree-level cross section leading to a correction of 
about 15% to the cross sections. The hght quark (effective) masses, are chosen such as to 
reproduce the SM. value of a''^{M^) — 127.77 including the light fermions contribution, 
which give the ~ 15% corrections compared to the use of a{0). For the SAi input masses 
see the last papers of Ref. [19] with the exception of rritop = 170.9 GeV. We will keep this 
rescaling in mind. This effect can be reabsorbed by using a scheme where the effective 
Q;(M|)eff. is used as input. 

ii) The Higgs sector: The Higgs sector is conceptually the trickiest. First we take Ma 
the pseudoscalar Higgs mass as an input parameter and require vanishing tadpoles. The 
extraction and definition of the ubiquitous tan/?, which at tree-level is identified as the 
ratio of the vev of the two Higgs doublet is the tricky part. Most schemes define the 
tan P counterterm at one-loop from a non physical quantity, such as the A^Z transition 
two-point function at — M\. It has become customary to take a DR definition, by 
only taking into account the "universal" ultraviolet part from such quantities, leaving 
out all finite parts. These prescriptions are however not gauge invariant, see for exam- 
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pie [34]. Moreover the "universal" part is only universal in the usual linear gauge. With 
the non-linear gauge fixing we implement, see section [273| our results would not be gauge 
invariant and one has to be very careful with the Ward identities. We leave this important 
issue to a forthcoming paper [35]. Nonetheless to conform with this widespread general 
usage, we also implement a DR scheme defined from a physical quantity, to be discussed 
shortly, which reproduces the usual counterterms defined from other quantities in the 
linear gauge. As known, the others Higgs masses m/^o (for the lightest CP-even), muo 
(for the heaviest CP-even) and m/^± (for the charged) receive corrections that can be very 
important. To be able to stick with an on-shell definition and in order to weigh the effect 
of the tan/5 scheme dependence, we also define two other schemes. One is based on the 
use of — s> T^T^ partial width to which the QED corrections have been extracted, we 
will refer to this scheme as the Arr scheme. For the third one, we take itiho as an input 
parameter and trade it for "tan /3" hence loosing one prediction, we will call this scheme 
MH. This scheme is also used in [29]. With tan (3 fixed, we can turn to the other sectors. 

iii) Neutralino and charginos: For the neutralino and chargino sector, we implement an 
on-shell scheme taking as input parameters the masses of the two charginos (this defines 
the counterterms to the SU(2) gaugino, wino mass M2 and to the higgsino, /i, parame- 
ter /i) and the mass of the LSP m^o (which completes the extraction of the U(l) gaugino, 

bino 6, mass Mi). The other neutralino masses m^o, m^o and m^o receive corrections 

' X2 ' A3 X4 

to their tree-level value. Obtaining finite corrections for the masses and decays is a not 
trivial test of the procedure. Here our implementation is quite similar to the one adopted 
in [36] when one sticks to the DR tan (3 . 

iv) Sfermions: For the slepton sector we use as input parameters the physical masses of 
the two charged sleptons which in the case of no-mixing define the i?-slepton soft breaking 
mass, and the SU{2) mass, Mj^, giving a correction to the sneutrino mass at one- 
loop. In the case of mixing one needs to fix the counterterm to the tri-linear coupling. 
The best option would have been to define this from a decay such as f2 fiZ. In the 
present letter we take a much simpler prescription, we solve the system by taking as input 
all three slepton masses. For the squark sector, for each generation three physical masses 
are needed as input to constrain the breaking parameters Mq^ for the SU{2) part, M^^, 
Mj^ for the i?-part. In case of mixing, the simplest prescription for the counterterms to 
the tri-linear couplings Au-, Ad derives from two conditions on the renormalised mixed 
two-point functions as is done in [37]. Our plan is to replace these conditions by an on- 
shell input such as a decay of the heavy squark to the lighter one and a Z, to conform 
with a fully on-shell scheme and study further the scheme dependence. 

Wave function renormalisation is introduced so that the residue at the pole of all physical 
particles is 1 and no-mixing is left between the different particles when on shell. This 
applies for all sectors. Dimensional reduction is used as implemented in the FFL bundle 
at one-loop through the equivalent constrained dimensional renormalisation [38]. Renor- 
malisation of the strong coupling constant and the gluino is not an issue for the examples 
we study here. 

We have verified our codes and schemes with different calculations on the market for 
a variety of correction to masses and other observables [28,34,36,37]. The code has 
also been used for corrections to SM. processes and also to one-loop induced processes 
XiXi ~^ 77; ^1:99 [11] relevant for indirect detection. 
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2.3 Non-Linear gauge- fixing 



We use a generalised non-linear gauge [11,39] adapted to the minimal supersymmetric 
model. The gauge fixing writes 

Cgf = -^\{d^-ie&-f^-igcwPZf,)W''+ + ^w^{v + 6h + uH + ikG3 + ipA)G^\^ 

^ -{d.Z + iz^{v + Ih + iE)Gzf - ^{p-lf- (1) 



Unlike the other parts of the model, Lqf is written in terms of renormalised fields and 
parameters. G3, are the Goldstone fields. We always work with ^■y,z,w = 1 so as to deal 
with the minimal set of loop tensor integrals. More details will be given elsewhere 



2.4 The diff'erent parts of the cross section 

The one-loop amplitudes consist of the virtual corrections ■4,^)^^^*^'^'^ and the counterterm 
contributions Act- ^iiYop^'^^ + -^ct niust be ultraviolet finite. On the other hand 
'^lioop^'^*"^ can contain infrared divergences due to photon and gluon virtual exchange. 
These are regulated by a small photon or gluon mass. For the QCD corrections we 
study here, this implementation does not pose a problem. The photon and gluon mass 
regulator contribution contained in the virtual correction should cancel exactly against 
the one present in the photon and gluon final state radiation. The photonic (gluonic) 
contribution is in fact split into a soft part, where the photon (gluon) energy is less than 
some small cut-off k^, A^°^^{E^^g < kc) and a hard part with A^'^g'''{E^^g > kc). The former 
requires a photon/gluon mass regulator. We use the usual universal factorised form with 
a simple rescaling for the case of the gluon correction. We take as = as{Mz) = 0.118. 



2.5 Checks on the calculation 

i) For each process and for each set of parameters, we first check the ultraviolet finiteness 
of the results. This test applies to the whole set of the virtual one-loop diagrams. The 
ultraviolet finiteness test is performed by varying the ultraviolet parameter Cuv = 1/f^- 
We vary Cuv by seven orders of magnitude with no change in the result. We content 
ourselves with double precision. 

ii) The test on the infrared finiteness is performed by including both the loop and the soft 
bremsstrahlung contributions and checking that there is no dependence on the fictitious 
photon mass A-^ or gluon mass A^. 

iii) Gauge parameter independence of the results is essential. It is performed through the 
set of the seven gauge fixing parameters defined in Eq. ([1]). The use of the seven parameters 
is not redundant as often these parameters check complementary sets of diagrams. It is 
important to note that in order to successfully achieve this test one should not include 
any width in the propagators. In fact our tree-level results do not include any width. 
Because of the parameters and the energies we consider, no width is required to regularise 
the cross sections. 
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iv) For the bremsstrahlung part we use VEGAS adaptive Monte Carlo integration package 
provided in the FFL bundle and verify the result of the cross section against CompHep [31]. 
We choose kc small enough and check the stability and independence of the result with 
respect to kc. 



2.6 Boltzmann equation, the small v expansion 

Having the collection of cross sections and the masses of the annihilating (and co-annihilatinj 
DM particles we could have passed them to micrOMEGAs for a very precise determination 
of the relic density based on a careful treatment of the Boltzmann equation. However, 
to weigh the impact of the corrections on the relic density it is worth to gain insight 
through an approximation in going from the cross sections to the relic density, especially 
that we have found these approximations to be, after all, rather excellent for the cases 
we study, including co-annihilations. Moreover corrections to the cross sections could 
be incorporated in the case of non-thermal production. All cross sections aij where i,j 
label the annihilating and co-annihilating DM particles are expanded in terms of the 
relative velocity Vij, which for neutralino annihilation is f = 2/3 = 2. /l — AuiIq/s. Away 



V ^i" 

from poles and thresholds, it is a very good approximation to write (TijVij = Qij + ^, 
keeping only the s-wave, ajj, and p-wave, hij coefficients. With T being the temperature, 
X = m^o/T, the thermal average gives 

{o-ij Vij) = aij + 6{bij - aij/A)/x, (2) 

with gi = 2 the neutralino spin degree of freedom (sdof), the co-annihilating particle of 
sdof Qi and mass contributes an effective relative weight of 

^i,e// = ^(1 + 5i)^^^exp(-x5i), 6i = {rrii - m^o) /m^o. (3) 

The total number of sdof is gejf = J2i9i,eff- good approximation for the relic density 
is obtained by carrying a simple one dimensional integration 

f 10 A 0.237 X 10-26cm3.s-i ^ P/ \ ^ / 2 

Me// = E^'^^^^^^--^^)- (4) 
ij 3eff 



Hj^bij that are needed to compute aij in Eq. (jlj) are given in cm^.s"^. xp represents 
the freeze-out temperature. g*{xp) is the effective degrees of freedom at freeze-out. 
is tabulated in micrOMEGAs and we read it from there. For the examples we will study 
a/ g*{xF) ~ 9.29. In the freeze-out approximation, xp can be solved iteratively from 



^21.2181 + In fe^^<"">'f^'l-'^. f"" c(c+2) . (5) 
1^ 10-26 lOOy j,(i-f)i^f 7 

where the neutralino mass is expressed in GeV. The numerical solutions of the density 
equation and hence the freeze-out suggest [6,40] c = 1.5 is a very good choice in most, but 
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not all, cases. Though we have verified that Eq. ([5]) converges quickly and agrees well with 
the result of micrOMEGAs, the results we give use Xp as extracted from micrOMEGAs. The 
loop corrected cross sections should also impact on the value of which is not exactly 
the same as the value extracted from the tree-level cross sections. However, the shift is 
marginal, though ultimately in a full computation this should be taken into account. Our 
results therefore use the same value oi xp at both tree and one-loop level. On the other 
hand to derive the relic density we rely, in this letter, on Eq. For The case of XiXi 
annihilation, the latter simplifies to 

- ( '^■^^'^ ^ lO^^^cm^.s-i 

V ^ / a + 3(6-a/4)/x^ 

The weight of a channel (see the percentages we will refer to later) corresponds to its 
relative contribution to J. 



2.7 Choosing points in the MSSM parameter space 

Current limits on the relic density, from WMAP and SDSS [3] give the 2a range 

0.092 < Qh^ < 0.124. (7) 

In this first exploratory study we thought it is best to consider different scenarii without 
worrying too much about the absolute value of the derived relic density in order to grasp 
the origin of the large corrections, if any. Our choice of scenarii was motivated by the 
physics issues, although our choice is biased towards the popular scenarii that emerge 
in mSUGRA within thermal production. Nonetheless our choice covers annihilation into 
fermions, gauge bosons and co-annihilation. This said, apart from the annihilation into 
gauge bosons, the derived relic density is either within this range or not overly outside. 
For the gauge bosons the motivation was to take a model that singles out the WW and 
ZZ final states channels. Moreover since the impact of the radiative corrections can be 
large there is not much sense in picking up a model based on its agreement with the 
current value of the relic density on the basis of a tree-level calculation. This said we have 
used micrOMEGAs as a guide, being careful about translations of effective couplings and 
input parameters. micrOMEGAs was also quite useful in justifying the approximations we 
use for deriving the relic density from the cross sections. We should also add that in this 
letter we do not apply the radiative corrections to all the subprocesses that can contribute 
to the relic density but only to those channels and subprocesses that contribute more than 
5% to the relic density. When calculating the correction to the relic density we include 
these channels at tree-level. 



3 Annihilation of a bino LSP into charged leptons 

The first example we take corresponds to the so-called bulk region, with a neutralino LSP 
which is mostly bino. The latter will couple with the particles of largest hypercharge, 
the R-sleptons. Therefore annihilation is into charged leptons. Because of the Majorana 
nature of the LSP, there is no s-wave contribution in the case of massless fermions. In 
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~0~0 lr>aW\ 

XiXi (36%) 


Tree 




A 


Dri 


1\ /f T T 


a 


0.081 




+38% 


+35% 


+15% 


b 


3.858 




+18% 


+ 18% 


+18% 




0.166 




0.138 


0.138 


0.141 








-17% 


-17% 


-15% 



Table 1: Bino case: Tree-level values of the s-wave (a) and p-wave (h) coefficients in 
units 10~^^cm^s~^, as well as the relative one-loop corrections in the DR, A^-j. and MH 
scheme. The derived relic density and its percentage change are also given. f36%j to the 
annihilation process refers to the percentage weight, at tree-level, of the channel. 



our case the contribution to the s-wave (at tree-level) is from the r's. In the radiation 
dominated standard scenario, to be consistent with the present value of the relic density 
, we take right sleptons as light as possible (but within the LEP limits) while all other 
particles (squarks, charginos, other neutralinos) heavy. The example we take has tan (3 = 
5. The relevant physical parameters are m^^o = 90.72 GeV, = 115.15 GeV and rrig^ = 
rrijx^ = 117.5 GeV. The masses of the charginos are m-+ = 200.64, m-+ = 610.47 GeV. To 
give an idea this reconstructs the tree-/ewe/ values for the gaugino and fiiggsino parameters 
as Ml = 90 GeV while M2 = 200, /i = —600 GeV leading to a neutralino which is almost 
100% bino: Xi = 0.9986 + O.OUw - O.OGShi + O.OOSh- 

The other physical input masses reconstruct the tree-level parameters rriej^ = llOGeV, 
M^^ = 250 GeV, M3 = Mq^ = m- = 800 GeV, Ma = 500 GeV and the tri-hnear 
coupling Af = 0. The contribution to the relic density is, then as expected, into leptons 
(98%) with the proportions as shown in Tabled! The e"'"e~, /i"*"/!" channels contribute for 
31% each. The difference between the three channels is accounted for by the contribution 
of the s-wave of the r final states and very little from the fact that fi is slightly lighter 
than the lightest smuon and selectron. 

Let us first comment on the p-wave contribution which gives the bulk of the contribu- 
tion to the relic density. The total correction is about 18% in this case. It is tempting to 
parameterize the corrections. In fact, had we used the value of the gauge coupling not at 
low energy but at the scale of the x? mass of order Mz the bulk of the correction would be 
absorbed. Indeed, (137.04/127.77)^ — 1 ~ 15%. In the few other examples we have looked 
at concerning the annihilation into leptons, we arrive at the same order of correction, see 
for example the corrections to the p-wave in the case of co-annihilation and even where 
there is some higgsino component as in the case of annihilation into bb. The other common 
trend is that the correction does not show any dependence on the tan (3 scheme we choose 
when there is a large bino component. This is not the case for the s-wave. In particular 
the MH scheme differs from the DR and Arr- The tan/3 dependence comes essentially 
from the Yukawa contribution, see Eq. (4.8) of ref [1]. The latter is also sensitive to the 
higgsino component of the neutralino that is also affected by the tan/3 change. The 
effect is more obvious in the case of scenario iii) , see Table |H Note however that even in 
the case of massless fermions, there is (though small) a contribution to the s-wave due to 
hard photon radiation. Hard photon radiation in association with a light charged fermion 
pair is not subject, for the s-wave amplitude, to the known helicity suppression when no 
photon is emitted [12]. Taking both the s and p-wave contribution leads to a correction 
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on the relic density of about 17%. As discussed a few lines above, using aeff{Mz) reduces 
the correction to the level of a few percent. 

4 Neutralino f co-annihilation 

In this scenario the LSP is still the lightest neutralino and we take it to be essentially 
bino though with a small higgsino component with a composition Xi = 0.9866 — 0.049ty + 
0.144/ii — 0.070/^2 and mass m^o = 162.34 GeV. We consider a scenario where the NLSP 
is the lightest fi with mass m^^ = 168.42 GeV coming mainly from its component. 
The lightest smuon and selectron are given masses so that they are thermodynamically 
irrelevant, they have a mass of about 195 GeV. Apart from = ITeV, for the squark 
and Higgs sector, as well as tan/3 , the parameters are the same as in the example in 
the bulk region. We want therefore to concentrate on co-annihilation involving only fi. 
With a mass difference between the LSP and NLSP of 6.08 GeV which is about only 
4% of the mass of the LSP, Xi^i annihilation is quite efficient with the two channels 
Xiff^ — > T^7, T^Z accounting for as much as half of the total contribution, see Table |2l 
f^f^ takes up a quarter of the total. Neutralino annihilation makes up for about 15%. 
The rest is made up by f^f^ — > 77, 'yZ. It is interesting to note that f^f^ ^ f^f+. Our 
approximation based on Eq. (jll) gives Qh'^ = 0.128, with xp = 26.5. 
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DR 


MH 
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0.002 




+200% 


+200% 


+200% 
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1.717 




+18% 


+19% 


+ 19% 
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MH 
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4.342 




+9% 


+9% 


+9% 


b 


-1.116 




+9% 


+8% 


+9% 


X^rf ^r±Z(10%) 


Tree 




Att 


DR 


MH 
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1.093 




+21% 


+21% 


+21% 


b 


-0.214 




+19% 


+19% 


+ 19% 


f±f± ^ r^r^ (23%) 


Tree 




Att 


DR 


MH 
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43.345 




+17% 


+ 17% 


+ 17% 


b 


-14.445 




+13% 


+13% 


+ 14% 


c 







-0.994 


-0.994 


-0.994 


Qh^ 


0.128 




0.117 


0.117 


0.117 


snh^ 






-9% 


-9% 


-9% 



Table 2: Neutralino Ti co- annihilation: meaning of the different cells is as in Table [H 
Percentages in the first column represent the weight of the corresponding channels. 

Compared to the bino-case in the bulk region where this channel accounted for the totality 
of the relic density, here it only makes up for 6%. Nonetheless the effect of radiative 
corrections on this channel are very similar to what we found in the scenario of section [31 
One may be misled by interpreting the 200% relative correction as a large correction to 
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the relic density. This is a relatively large correction to the s-wave contribution, but in 
absolute terms this correction is totally negligible compared to the correction brought 
about by the p-wave contribution at the cross section level as well as after taking the 
thermal average, notwithstanding that the whole channel in the co-annihilation region 
has little weight which further dilutes the effect of such large relative correction. As 
pointed out in the previous section this 200% relative correction is due to the smallness 
of the s-wave XiXi ~* ff which is offset by the hard photon emission that now allows 
for an s-wave contribution [12]. As discussed in the first reference in [12], the relative 
importance of hard radiation reduces fast as the mass of the intermediate slepton increases. 
This explains why the relative effect is more prominent in the co-annihilation region. 

Before getting into the details about each of the main contributing channels that 
involve co-annihilation, let us point at some of their common features. The bulk of 
the contribution comes now from the s-wave, especially after taking into account thermal 
averaging. Another common feature is that the tan/? scheme dependence is hardly 
noticeable. Moreover, the corrections to the s- and p-wave are, within a margin of 4%, 
the same. 

For Xi'^^i^ ~^ '^^1 the correction to the p-wave is within only a 1% from the correction to 
the s-wave. Here the electroweak correction amounts to about 9%. This is about half the 
correction we find in all other channels. The reason is the following, the effective coupling 
for the emitted photon should still be taken at = 0, and therefore effectively since there 
is only one neutralino coupling, this should be taken at the scale of rn^o. Proper use of 
the effective couplings here absorbs about 7% of the correction, leaving therefore a 2% 
correction. For the tZ final state the use of the effective coupling would leave out a 6% 
correction to the s-wave after absorbing the correction due to the effective coupling. 



The radiative correction to f^f^ — >■ r^r^ reveals a quite interesting feature as can be 
seen from Fig. [1] which shows the dependence of the cross section as a function of the 
relative velocity (or rather the square of it). It is from this velocity dependence that we 
usually extract the values of the coefficients of the s- and p-wave contribution, at both 
the tree-level and one- loop, that we need to calculate the relic density. The figure extends 
to v'^ values as large as 0.5 while we could have contented ourselves with a maximal value 
of = 0.3, considering the typical value that one obtains for the freeze-out temperature 
{v"^) ^ 6/xf ^ 0.2. Even so, the figure shows that in this case a fit to the tree-level cross 
section in the form a + hv"^ works quite well. For the one-loop correction a polynomial fit 
does not do for low enough velocities. There is a large negative correction for f ^ 0. This 
correction is in fact very easy to understand. It is the perturbative one-loop manifestation 
of the non relativistic Coulomb-Sommerfeld effect [41]. With the tree-level cross section 
denoted as ctq and a^v = oq + hov"^, the one-loop perturbative cross section for the same- 
sign stau annihilation f^f^ , '^coui°''' ^Vidi that 

1— loop 

^Coul. _ (g) 
(To V 

We thus expect for the one-loop cross section, cxi 

aiv = ai + hiv'^ + Ci/v with ci = — vraao. (9) 
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30=43.345 
bo=-1 4.445 
3^=50.803 
bi =-16.377 
Ci =-0.994 
32=50.789 
=-16.286 
:-0.991 



0.1 
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One Loop 

0.4 



0.5 



Figure 1: f^f^ r r as a function of the square of the relative velocity both for tree- 
level (solid line) and at one-loop (dotted line) in the A^r scheme. Fits to the s-wave (a), 
p-wave (b) and the Coulomb factor are given. ao,6o for tree-level, ai,bi for the loop with 
Ci given by Eq. fli, 6i, C2 when all three parameters are extracted from the fit. The two 
fits at one-loop are not distinguishable in the figure. 



Fig. [T] reflects this repulsive behaviour perfectly. In fact we made a fit to the one-loop 
result with the function aiv = ai + hiv'^ + ci/f , first with ci as given in Eq. ([9]) and then 
with Ci not constrained. The two fits are practically indistinguishable in Fig. [TJ Our 
automatic calculation code captures this effect perfectly. 

One may ask about how to deal with the 1/v singularity. In fact when calculating the 
relic density, the one- loop 1/v singularity at the level of the cross section is tamed after 
thermal average, oc j^{(^v 1^^) (erf), see also Eq. [21 At the end its contribution 

to the relic density compared to a\ is approximately x^{a\ — 2ai^(x^Txxp). In words, non 
zero temperature of the problem provides a cut-off. One can also ask whether the one-loop 
result from the Coulomb- Sommerfeld effect is sufficient. As seen, the QED correction is 
of (9(7ra/f) ~ 0.17 with v typical of the freeze-out temperature. Therefore in our case a 
one-loop treatment seems to be sufficient especially that the t^t^ r^r^ is not the most 
dominant channel. This said, these non-relativistic QED \jv threshold corrections can 
be resummed to all orders. This resummation as originally performed by Sommerfeld [41] 
has been known for quite a long time in quantum mechanics, see [42] for a textbook 
treatment, and amounts to solving the Schrodinger equation in the Coulomb potential. 
With Xnr = —2TTa/v for the same-sign stau annihilation, the s-wave factor resums to 

S„, = ^^ = l-- + i(-)V- (10) 

1 — e "^'"^ V \ V / 
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One might question the vahdity of Eq. [TO] in our case where fi is not stable. Finite 
decay width can of course act as a cut-off for the l/t> corrections, see the case of W pair 
production [43] or slepton pair production at threshold [44]. In our case width effects 
are of no importance since the characteristic time of the Coulomb interaction, 
typical of velocities at freeze-out is much smaller than the decay time, l/F^^, since in our 
example Ffj/mfj ^ 2 10^^ and v"^ ~ 0.2. For smaller 6m = m^^ — m^o, the width effect is 
even more negligible whereas for larger 6m, the fifi channel would be thermodynamically 
irrelevant. Therefore in our particular case the resummation can be taken from the old 
Sommerfeld result. Nonetheless, especially after thermal averaging, in our case this type 
of QED correction is well approximated by the one-loop approximatioiil. 

Taking now all the effects and contributions in our specific example we find an overall 
correction of —8.6% to the relic density with a corrected value of Qh'^ = 0.117. As we 
can see this value would not have been approached with a naive overall rescaling of the 
effective couplings. Nonetheless, apart from some 4% correction, most of the effect seems 
to be explained in terms of a proper usage of effective couplings and the Coulomb effect. 
In fact in the total contribution, the Coulomb effect is diluted and changes the results for 
the relic density by about only 1.5%. 

5 Annihilation of a mixed gaugino-higgsino LSP into 
vector bosons 

Having studied annihilation into fermions, annihilation into the weak vector bosons is 
quite interesting. In the context of mSUGRA this occurs for example in the so-called 
focus point region. In order not to mix issues, we do not consider in this letter a scenario 
where the LSP neutralino is either dominantly higgsino or wino, therefore avoiding that 
XiXt annihilation is of relevance in this case. We seek a scenario with a neutralino where 
the largest component is still bino but where one has a substantial higgsino and wino 
component. In our example one has Xi = 0.8196 — 0.231w — 0.470/ii — 0.232/i2 with 
m^o = 102.89, m-+ = 125.13, m-+ = 215.27 GeV. All other masses (outside the chargino- 
neutralino sector) are taken to be very heavy at ITeV. This is also to avoid contamination 
from annihilation into fermions. It would however be worth to study the impact of the 
sfermions on the radiative corrections. In this case we have not made the extra effort of 
searching for a set with a relic density within the WMAP range. Here annihilation into 
WW and ZZ accounts for 80%, see Table [3] with a few other channels below 2% each. 
These involve XiXt i^ito light quarks or WZ which we take into account only at tree-level. 

First of all, we see that the corrections which affect annihilation into ZZ and WW 
are about the same (within 5% in the 3 tan/5 schemes). Moreover the correction to the 

^ In a situation, which is not the case here, where the r^f^ and f^f ^ would contribute equally at tree- 
level, the Coulomb- Sommerfeld correction would cancel after adding the two-channels. The correction in 
f^f^ would be attractive and given by changing the sign of Xnr in Eg. llOI The non-perturbative effects 
of the Coulomb-Sommerfeld-like corrections that might occur when coloured states are involved [9,45] 
need more care because of the strong QCD coupling, bound states effects might be relevant. The effect of 
the latter is even more important with models with TeV and multi TeV dark matter almost degenerate 
with a charged component [16,46]. The non pcrturbativc effects on indirect detection in these models is 
even more dramatic [15,46]. 
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Tree 






DR 


MH 


a 


3.099 




-27% 


-2% 


+44% 


D 


o.yoi 




-oz/o 


- r /o 


+ 00/0 


XiXi ^ (5%) 


Tree 




/I 


UK 


MH 


a 


0.159 




-22% 


+3% 


+50% 


b 


0.787 




-30% 


-6% 


+39% 




0.053 




0.068 


0.054 


0.039 


snh'-' 






+28% 


+2% 


-26% 



Table 3: Mixed case: As in TableUi Note again that the percentages in the first column 
next to the process refer to the percentage weight, at tree-level, of that particular channel. 



s-wave and p-wave are of the same order, see also Fig. [2] for the v dependence of the 
ZZ cross section and the extraction of the s- and p-wave coefficients. However this is 
not the most important conclusion. The most important lesson is that there is a very 
large tan/3 scheme dependence. In some other investigation concerning Higgs decays, 
we had noticed, as was also pointed in [34] with a similar scheme, that the MH scheme 
can lead to very large corrections. However in many instances, like what we saw in the 
case of the bino annihilation or co-annihilation, the 5tan/5 is screened. Unfortunately 
in this mixed neutralino scenario, the 6 tan (3 dependence can be potentially enhanced by 
1 / (/i^ — M|) from the renormalisation of for example. This needs further investigation. 
In this model, already at tree-level, we had noticed that the cross sections were very 
sensitive to small changes in the underlying parameters. Apart from the tan/3 scheme 
dependence, the corrections in the DR scheme look modest especially for the dominant s- 
wave contribution. However one should not forget that these small corrections are within 
the use of a in the Thomson limit. Switching to a scale of order Mz, the corrections 
are large of order 20%. Therefore our conclusion is that in such a scenario there are 
genuine large corrections in all three schemes we have considered. This also confirms the 
study of the chargino neutralino system at one-loop in [47] which showed that though the 
corrections to the masses are modest, there can be a large (of order 30%) change in the 
gaugino-higgsino component and hence a large impact on cross sections. 



6 Annihilation of neutralinos to bb and a not too 
heavy pseudo-scalar Higgs 

We expose this last example for illustrative purposes. Indeed the one-loop perturbative 
treatment of the Higgs coupling to a bottom pair using the bottom pole mass, here we 
have taken all along = 4.7 GeV, is far from describing the bulk of the radiative QCD 
corrections which as we know need to be resummed both for the running effective 
purely from QCD and from the so-called Amb effects. The latter being more impor- 
tant for high tan/? . These effects are already taken in micrOMEGAs [6] for example. 
The purpose here therefore is to see whether there are other possible effects, though 
smaller, that are captured in a complete one-loop calculation. Ideally one would like 
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Figure 2: XiXi ~^ as a function of the square of the relative velocity both for tree-level 
(solid line) and at one-loop (dotted line) in the DR scheme. Fits to the s-wave (a), p-wave 
(b) are indicated. 



to subtract the known universal one-loop QCD corrections from the full one-loop QCD 
corrections that can occur for example from box diagrams, for v ^ 0, outside the Higgs 
resonance. By the way because of the Majorana nature of the LSP, at the smallest ve- 
locities the most important Higgs resonance is the pseudo-scalar Higgs, A. In any case 
for a precise calculation of the relic density, non resonant contributions should be taken 
into account as thermal average is to applied and would bring some smearing. There- 
fore, here we concentrate on XiXi ~^ bb where the scalar resonance is not negligible. 
Again like we argued for the previous example we have to rely here also on some higgsino 
component. The composition of the LSP is Xi = 0.8966 — 0.161w — 0.378/;.i — O.I66/12 
with m^o = 105.74, m-+ = 130.99, = 225.10 GeV. At tree-level the system is 
(re) constructed from a set Mi,M2,/i = 110,150,-180 GeV and tan/3 = 5. Compared 
to the previous case where all other masses were around ITeV, to bring out the effect of 
the 6's in the final state we lower first Ma such that Ma = 300 GeV. The masses of all 
sfermions are around 200 GeV for the dominant fn and 250 GeV for the dominant f^. 
The mass of the gluino is also lowered to be 400 GeV. 

At tree- level with rrif, = 4.7 GeV, the dominant modes are annihilation into WW followed 
by bb (about 10% smaller), rr and ZZ are about the same level but a factor 10 smaller 
than the WW channel. We show, see Table HI the WW and ZZ channel in order to make a 
comparison with the previous case of a mixed bino with a substantial higgsino component. 
Here the tan /? scheme dependence has considerably reduced especially between the Arr 
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MH 
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-1% 


+3% 
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5a/a QCD 


-26% 


-26% 


-26% 


56/6 EW 


-1% 


+3% 


+29% 


db/b QCD 


-30% 


-30% 


-30% 



Table 4: Mixed case 2: As in TablelM for the array on the left. The array on the right 
gives the relative corrections to the bb channel for the QCD and EW corrections. 



and DR scheme. Note also that in the case of rr channel there is also a discrepancy with 
the MH scheme and the other two schemes compared to the almost pure bino case. This 
again is due to the larger contribution of the higgsino tan [3 dependent part, naturally in 
the A-exchange not present in the pure bino case but also in the f exchange. 

Let us look at what we obtain for the hb channel. The electroweak corrections do show 
some tan/3 scheme dependence. Compared to the electroweak corrections in the DR 
scheme of tan (3 , the QCD corrections are larger by an order of magnitude, they amount 
to about —30% in both the s- and p-wave. If one assumes that most of these corrections 
arise from the A ^ hb vertex, then we know that there are large logs resulting from the 
anomalous dimension of the pseudoscalar (and for that matter the scalar). Its one loop 
part can be found in [48] and amounts to 




There is also an important SUSY QCD correction termed Amb [49], see section B.4.2 in 
the second paper of [6]. In this scenario, adding these two corrections amounts to about 
—35%. Subtracting these from the correction we calculate for the s-wave, leaves us with 
about +9% QCD corrections. The QCD corrections from the anomalous dimension and 
the Amb effect were extracted for the known effect to A — > hb. However, since at f = 0, 
for the s-wave cross section, the neutralino system constructs a pseudoscalar because 
of its Majorana nature, the same corrections should affect the a coefficient even if the 
contribution from the pseudo-scalar Higgs is negligible. 



7 Conclusions 



We have performed the first electroweak corrections to some important processes relevant 
for the relic density of neutralinos in supersymmetry. This has become possible thanks 



15 



to an automated code for the calculation of loop corrections in the MSSM that will allow 
to perform with the same tools and the same framework (scheme dependence,..) analyses 
at one-loop at the collider and for dark matter. Our findings suggest that in the case 
of a dominantly bino neutralino, a large part of the correction can be accounted for 
through an effective electromagnetic coupling at the scale of the neutralino mass. Even 
so, complete one loop corrections would be needed to match the foreseen precision of 
PLANCK. The corrections to the relic density are not sensitive much to the way tan/3 is 
renormalised. In the case of co-annihilation of a bino and stau, the conclusion is similar 
but one has to be wary of possible Coulomb-Sommerfeld corrections. For a neutralino 
LSP which is strongly mixed, the corrections are large and the tan (5 scheme dependence 
not neghgible at all. More investigation of such scenarii should be conducted. Some 
QCD (and SUSY QCD) corrections affecting final states quarks in the case of neutralino 
annihilation need that one goes beyond one-loop. Some of these corrections have been 
identified and already implemented in a code such as micrOMEGAs. Apart from these 
corrections, there remain some additional one-loop corrections that should be taken into 
account. Before generalising these conclusions, more work is needed. However, the tools 
exist. The next step is to interface our code for the loop calculations with a dedicated 
relic density calculator, avoiding double counting of some of the one-loop corrections 
implemented as effective operators in the relic density calculator. Work in this direction 
has already begun based on micrOMEGAs as the relic density calculator. 
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